analysis_tag <- '01_jaeger_descriptive'
library(pheatmap)
library(knitr)
library(scater)
library(Seurat)
library(Cairo)
library(biomaRt)
library(dplyr)
## library(reshape2)
#library(iSEE)
library(Matrix)
library(scran)
## library(mvoutlier)
## library(shiny)
## library(kableExtra)
## library(velocyto.R)
## library(cowplot)
## library(corrplot)
library(DT)
library(viridis)
## library("CellMixS")
ac <- function(col, alpha = 1) {
apply(sapply(col, col2rgb)/255, 2, function(x) rgb(x[1],
x[2], x[3], alpha = alpha))
}
Raw counts retrieval. Assumes all matrices share gene id and gene names
fns <- list.files(file.path("..", "data"), "*gz")
ec <- list()
for (fn in fns) {
tag <- gsub("_rawCount.txt.gz", "", basename(fn))
ec[[tag]] <- read.table(file.path("..", "data", fn), sep = "\t",
header = TRUE)
rownames(ec[[tag]]) <- sprintf("%s_%s", ec[[tag]]$gene_id,
ec[[tag]]$gene_name)
gene_ids <- ec[[tag]]$gene_id
ec[[tag]] <- ec[[tag]][, -c(1, 2)]
}
sapply(ec, function(x) head(colnames(x)))
Ascl1_A_12wk Ascl1_A_5d Ascl1_B_12wk
[1,] "Ascl1_12wk_A01" "Ascl1_d4_plate1_A01" "Ascl1_12wk_2_20190514BJ_A01"
[2,] "Ascl1_12wk_A02" "Ascl1_d4_plate1_A02" "Ascl1_12wk_2_20190514BJ_A02"
[3,] "Ascl1_12wk_A03" "Ascl1_d4_plate1_A03" "Ascl1_12wk_2_20190514BJ_A03"
[4,] "Ascl1_12wk_A04" "Ascl1_d4_plate1_A04" "Ascl1_12wk_2_20190514BJ_A04"
[5,] "Ascl1_12wk_A05" "Ascl1_d4_plate1_A05" "Ascl1_12wk_2_20190514BJ_A05"
[6,] "Ascl1_12wk_A06" "Ascl1_d4_plate1_A06" "Ascl1_12wk_2_20190514BJ_A06"
Ascl1_B_5d Gli1_A_12wk
[1,] "Ascl1_5d_II_20190403BJ_A01" "Gli1_12wk_20180822BJ_A01"
[2,] "Ascl1_5d_II_20190403BJ_A02" "Gli1_12wk_20180822BJ_A02"
[3,] "Ascl1_5d_II_20190403BJ_A03" "Gli1_12wk_20180822BJ_A03"
[4,] "Ascl1_5d_II_20190403BJ_A04" "Gli1_12wk_20180822BJ_A04"
[5,] "Ascl1_5d_II_20190403BJ_A05" "Gli1_12wk_20180822BJ_A05"
[6,] "Ascl1_5d_II_20190403BJ_A06" "Gli1_12wk_20180822BJ_A06"
Gli1_A_5d Gli1_B_5d
[1,] "Gli1_20180503BJ_Plate_I_B13" "X20190219_Gli1_5d_A01"
[2,] "Gli1_20180503BJ_Plate_I_B14" "X20190219_Gli1_5d_A02"
[3,] "Gli1_20180503BJ_Plate_I_B15" "X20190219_Gli1_5d_A03"
[4,] "Gli1_20180503BJ_Plate_I_B16" "X20190219_Gli1_5d_A04"
[5,] "Gli1_20180503BJ_Plate_I_B18" "X20190219_Gli1_5d_A05"
[6,] "Gli1_20180503BJ_Plate_I_B20" "X20190219_Gli1_5d_A06"
Checking all rownames are equivalent
for (item in names(ec)) {
for (item2 in names(ec)) {
stopifnot(all(rownames(ec[[item]]) == rownames(ec[[item2]])))
}
}
Metadata (currently just the plate)
meta <- list()
for (tag in names(ec)) {
meta[[tag]] <- data.frame(cell = colnames(ec[[tag]]), plate = tag)
rownames(meta[[tag]]) <- meta[[tag]]$cell
}
sce <- list()
for (item in names(ec)) {
stopifnot(all(rownames(meta[[item]]) == colnames(ec[[item]])))
sce[[item]] <- SingleCellExperiment(list(counts = as.matrix(ec[[item]])),
colData = meta[[item]])
}
Checking all rownames are identical
for (item in names(sce)) {
for (item2 in names(sce)) {
stopifnot(all(rownames(sce[[item]]) == rownames(sce[[item2]])))
}
}
Getting gene feature annotation. Assumes all data matrices share the same ensembl GIDs, and that this is mouse data. Only queries it the first time.
if (!file.exists(file.path("..", "data", "gene_annotation.RData"))) {
## ids = gene_ids
filters = "ensembl_gene_id"
attributes = c(filters, "chromosome_name", "gene_biotype",
"start_position", "end_position")
biomart = "ENSEMBL_MART_ENSEMBL"
dataset = "mmusculus_gene_ensembl"
host = "www.ensembl.org"
bmart <- biomaRt::useMart(biomart = biomart, dataset = dataset,
host = host)
feature_info <- biomaRt::getBM(attributes = attributes, filters = filters,
values = gene_ids, mart = bmart)
mm <- match(gene_ids, feature_info[[filters]])
feature_info_full <- feature_info[mm, ]
save(feature_info_full, file = file.path("..", "data", "gene_annotation.RData"))
} else {
load(file.path("..", "data", "gene_annotation.RData"))
}
for (item in names(sce)) {
object <- sce[[item]]
old_rdata <- rowData(object)
keep_cols <- !(colnames(old_rdata) %in% colnames(feature_info_full))
new_rdata <- cbind(old_rdata[, keep_cols], feature_info_full)
rowData(object) <- new_rdata
rowData(object)$ensembl_gene_id <- gene_ids
sce[[item]] <- object
rm(object, old_rdata, keep_cols, new_rdata)
}
Looking for mitochondrial/spike ins for filtering out cells afterwards.
We exclude genes that are not expressed in any cell.
for (item in names(sce)) {
sce[[item]] <- sce[[item]][rowSums(counts(sce[[item]]) >
0) > 0, ]
is.mt <- grepl("MT", rowData(sce[[item]])$chromosome_name)
is.spike <- grepl("ERCC-", rownames(rowData(sce[[item]])))
isSpike(sce[[item]], "Spike") <- is.spike
dim(sce[[item]])
## sce[[item]] <- normalize(sce[[item]])
sce[[item]] <- calculateQCMetrics(sce[[item]], feature_controls = list(ERCC = is.spike,
mito = is.mt))
}
Mind the bump in library sizes, are these cells duplets or different anyhow?
for (item in names(sce)) {
print(item)
par(mfrow = c(2, 2), mar = c(5.1, 4.1, 0.1, 0.1))
hist(sce[[item]]$total_counts/1000, xlab = "Library sizes (thousands)",
main = "", breaks = 20, col = "grey80", ylab = "Number of cells")
hist(sce[[item]]$total_features_by_counts, xlab = "Number of expressed genes",
main = "", breaks = 20, col = "grey80", ylab = "Number of cells")
hist(sce[[item]]$pct_counts_Spike, xlab = "ERCC proportion (%)",
ylab = "Number of cells", breaks = 20, main = "", col = "grey80")
hist(sce[[item]]$pct_counts_mito, xlab = "Mitochondrial proportion (%)",
ylab = "Number of cells", breaks = 20, main = "", col = "grey80")
}
[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"
for (item in names(sce)) {
print(item)
par(mfrow = c(2, 2), mar = c(5.1, 4.1, 1.1, 1.1))
plot(density(sce[[item]]$total_counts/1000), xlab = "library sizes (thousands)",
main = "")
rug(sce[[item]]$total_counts/1000)
plot(y = sce[[item]]$total_counts/1000, x = sce[[item]]$pct_counts_mito,
pch = 20, ylab = "library size (thousands)", xlab = "mitochondrial proportion (%)",
col = ac("black", 0.5))
plot(y = sce[[item]]$total_counts/1000, x = sce[[item]]$total_features_by_counts,
pch = 20, ylab = "library size (thousands)", xlab = "number of genes",
col = ac("black", 0.5))
plot(y = sce[[item]]$pct_counts_mito, x = sce[[item]]$total_features_by_counts,
pch = 20, xlab = "number of genes", ylab = "mitochondrial proportion (%)",
col = ac("black", 0.5))
}
[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"
Removal of outliers for the library size and the number of expressed features and mitochondrial proportions.
for (item in names(sce)) {
print(item)
libsize.drop <- isOutlier(sce[[item]]$total_counts, nmads = 3,
type = "lower", log = TRUE)
## feature.drop <-
## isOutlier(sce[[item]]$total_features_by_counts, nmads=3,
## type='lower', log=TRUE)
feature.drop <- sce[[item]]$total_features_by_counts >= 5000 |
sce[[item]]$total_features_by_counts < 1500
spike.drop <- isOutlier(sce[[item]]$pct_counts_Spike, nmads = 3,
type = "higher")
mito.drop <- isOutlier(sce[[item]]$pct_counts_mito, nmads = 2,
type = "higher")
manualmito.drop <- sce[[item]]$pct_counts_mito >= 5
## extra libsize outlier (manual) manuallibsize.drop <-
## sce[[item]]$total_counts < 100000
sce[[item]] <- sce[[item]][, !(libsize.drop | feature.drop |
spike.drop | mito.drop | manualmito.drop)]
print(data.frame(ByLibSize = sum(libsize.drop), ByFeature = sum(feature.drop),
BySpike = sum(spike.drop), ByMito = sum(mito.drop), ByManualMitoDrop = sum(manualmito.drop[!mito.drop]),
Remaining = ncol(sce[[item]])))
}
[1] "Ascl1_A_12wk"
ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1 1 309 0 19 269 60
[1] "Ascl1_A_5d"
ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1 24 31 0 43 0 308
[1] "Ascl1_B_12wk"
ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1 49 71 0 24 0 291
[1] "Ascl1_B_5d"
ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1 62 27 0 11 0 297
[1] "Gli1_A_12wk"
ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1 50 171 0 37 0 195
[1] "Gli1_A_5d"
ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1 8 26 0 41 8 204
[1] "Gli1_B_5d"
ByLibSize ByFeature BySpike ByMito ByManualMitoDrop Remaining
1 44 215 0 67 54 155
Cell QC (number of genes detected, number of reads, mitochondrial) after outlier removal.
for (item in names(sce)) {
print(item)
par(mfrow = c(2, 2), mar = c(5.1, 4.1, 0.1, 0.1))
hist(sce[[item]]$total_counts/1000, xlab = "Library sizes (thousands)",
main = "", breaks = 20, col = "grey80", ylab = "Number of cells")
hist(sce[[item]]$total_features_by_counts, xlab = "Number of expressed genes",
main = "", breaks = 20, col = "grey80", ylab = "Number of cells")
hist(sce[[item]]$pct_counts_Spike, xlab = "ERCC proportion (%)",
ylab = "Number of cells", breaks = 20, main = "", col = "grey80")
hist(sce[[item]]$pct_counts_mito, xlab = "Mitochondrial proportion (%)",
ylab = "Number of cells", breaks = 20, main = "", col = "grey80")
}
[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"
for (item in names(sce)) {
print(item)
par(mfrow = c(2, 2), mar = c(5.1, 4.1, 1.1, 1.1))
plot(density(sce[[item]]$total_counts/1000), xlab = "library sizes (thousands)",
main = "")
rug(sce[[item]]$total_counts/1000)
plot(y = sce[[item]]$total_counts/1000, x = sce[[item]]$pct_counts_mito,
pch = 20, ylab = "library size (thousands)", xlab = "mitochondrial proportion (%)",
col = ac("black", 0.5))
plot(y = sce[[item]]$total_counts/1000, x = sce[[item]]$total_features_by_counts,
pch = 20, ylab = "library size (thousands)", xlab = "number of genes",
col = ac("black", 0.5))
plot(y = sce[[item]]$pct_counts_mito, x = sce[[item]]$total_features_by_counts,
pch = 20, xlab = "number of genes", ylab = "mitochondrial proportion (%)",
col = ac("black", 0.5))
}
[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"
In this setting, feature controls are mitochondrial genes.
for (item in names(sce)) {
print(item)
print(plotHighestExprs(sce[[item]], exprs_values = "counts"))
p1 <- plotColData(sce[[item]], x = "total_counts", y = "total_features_by_counts")
p2 <- plotColData(sce[[item]], x = "pct_counts_feature_control",
y = "total_features_by_counts")
p3 <- plotColData(sce[[item]], x = "pct_counts_feature_control",
y = "pct_counts_in_top_50_features")
p4 <- plotColData(sce[[item]], x = "total_counts", y = "pct_counts_in_top_500_features")
multiplot(p1, p2, cols = 2)
multiplot(p3, p4, cols = 2)
}
[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"
for (item in names(sce)) {
tryCatch({
print(item)
genes_keep <- rowSums(counts(sce[[item]]) > 1) >= 10
sce[[item]] <- sce[[item]][genes_keep, ]
sce[[item]] <- scran::computeSumFactors(sce[[item]])
## no spike ins here sce[[item]] <-
## scran::computeSpikeFactors(sce[[item]])
sce[[item]] <- normalize(sce[[item]])
}, error = function(x) print(x))
}
[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"
Cell cycle phase estimation (mouse)
for (item in names(sce)) {
print(item)
mouse.pairs <- readRDS(system.file("exdata", "mouse_cycle_markers.rds",
package = "scran"))
assignments <- cyclone(sce[[item]], mouse.pairs, gene.names = rowData(sce[[item]])$ensembl_gene_id)
print(table(assignments$phase))
colData(sce[[item]])$est_cell_cycle_phase <- assignments$phases
}
[1] "Ascl1_A_12wk"
G1 G2M S
65 50 8
[1] "Ascl1_A_5d"
G1 G2M S
276 27 6
[1] "Ascl1_B_12wk"
G1 G2M S
297 8 2
[1] "Ascl1_B_5d"
G1 G2M S
289 17 6
[1] "Gli1_A_12wk"
G1 G2M S
241 53 4
[1] "Gli1_A_5d"
G1 G2M S
201 11 7
[1] "Gli1_B_5d"
G1 G2M S
250 21 11
Transformation to Seurat v3 objects. Please note that during data scaling we regress out both number of UMIs (in this case this means library size, because the protocol has no UMIs on it) and the mitochondrial proportions.
RANGE_PCS <- 1:30
## https://satijalab.org/seurat/pancreas_integration_label_transfer.html
so <- list()
for (item in names(sce)) {
print(item)
so[[item]] <- as.Seurat(sce[[item]])
so[[item]]@meta.data$orig.ident <- so[[item]]@meta.data$plate
Idents(so[[item]]) <- "plate"
}
[1] "Ascl1_A_12wk"
[1] "Ascl1_A_5d"
[1] "Ascl1_B_12wk"
[1] "Ascl1_B_5d"
[1] "Gli1_A_12wk"
[1] "Gli1_A_5d"
[1] "Gli1_B_5d"
First evaluating the unintegrated data structure
## so <- set.all.ident(so, id = identity)
merged <- merge(x = so[[1]], y = so[2:length(so)])
merged <- NormalizeData(object = merged)
merged <- FindVariableFeatures(object = merged)
merged <- ScaleData(object = merged)
merged <- RunPCA(object = merged)
PC_ 1
Positive: ENSMUSG00000006333-Rps9, ENSMUSG00000003657-Calb2, ENSMUSG00000042436-Mfap4, ENSMUSG00000037984-Neurod6, ENSMUSG00000023927-Satb1, ENSMUSG00000038642-Ctss, ENSMUSG00000024754-Tmem2, ENSMUSG00000036905-C1qb, ENSMUSG00000030707-Coro1a, ENSMUSG00000020135-Apc2
ENSMUSG00000026360-Rgs2, ENSMUSG00000025400-Tac2, ENSMUSG00000036887-C1qa, ENSMUSG00000036896-C1qc, ENSMUSG00000031822-Gse1, ENSMUSG00000028581-Laptm5, ENSMUSG00000058715-Fcer1g, ENSMUSG00000024621-Csf1r, ENSMUSG00000029771-Irf5, ENSMUSG00000030579-Tyrobp
ENSMUSG00000021876-Rnase4, ENSMUSG00000052336-Cx3cr1, ENSMUSG00000025511-Tspan4, ENSMUSG00000022619-Mapk8ip2, ENSMUSG00000031821-Gins2, ENSMUSG00000052160-Pld4, ENSMUSG00000074622-Mafb, ENSMUSG00000021423-Ly86, ENSMUSG00000031530-Dusp4, ENSMUSG00000007021-Syngr3
Negative: ENSMUSG00000026424-Gpr37l1, ENSMUSG00000020591-Ntsr2, ENSMUSG00000045092-S1pr1, ENSMUSG00000050953-Gja1, ENSMUSG00000030605-Mfge8, ENSMUSG00000006205-Htra1, ENSMUSG00000007097-Atp1a2, ENSMUSG00000022132-Cldn10, ENSMUSG00000023913-Pla2g7, ENSMUSG00000028517-Plpp3
ENSMUSG00000017390-Aldoc, ENSMUSG00000054252-Fgfr3, ENSMUSG00000029309-Sparcl1, ENSMUSG00000022037-Clu, ENSMUSG00000004892-Bcan, ENSMUSG00000029838-Ptn, ENSMUSG00000040055-Gjb6, ENSMUSG00000035805-Mlc1, ENSMUSG00000030428-Ttyh1, ENSMUSG00000041329-Atp1b2
ENSMUSG00000030495-Slc7a10, ENSMUSG00000055254-Ntrk2, ENSMUSG00000017009-Sdc4, ENSMUSG00000028128-F3, ENSMUSG00000024411-Aqp4, ENSMUSG00000032281-Acsbg1, ENSMUSG00000058135-Gstm1, ENSMUSG00000046240-Hepacam, ENSMUSG00000006522-Itih3, ENSMUSG00000005360-Slc1a3
PC_ 2
Positive: ENSMUSG00000036896-C1qc, ENSMUSG00000036887-C1qa, ENSMUSG00000036905-C1qb, ENSMUSG00000028581-Laptm5, ENSMUSG00000052336-Cx3cr1, ENSMUSG00000030579-Tyrobp, ENSMUSG00000038642-Ctss, ENSMUSG00000054675-Tmem119, ENSMUSG00000021665-Hexb, ENSMUSG00000015852-Fcrls
ENSMUSG00000023992-Trem2, ENSMUSG00000048163-Selplg, ENSMUSG00000058715-Fcer1g, ENSMUSG00000059498-Fcgr3, ENSMUSG00000024621-Csf1r, ENSMUSG00000036908-Unc93b1, ENSMUSG00000036353-P2ry12, ENSMUSG00000020101-Vsir, ENSMUSG00000051504-Siglech, ENSMUSG00000040747-Cd53
ENSMUSG00000027848-Olfml3, ENSMUSG00000021423-Ly86, ENSMUSG00000052160-Pld4, ENSMUSG00000018774-Cd68, ENSMUSG00000002111-Spi1, ENSMUSG00000029919-Hpgds, ENSMUSG00000036362-P2ry13, ENSMUSG00000021666-Gfm2, ENSMUSG00000015396-Cd83, ENSMUSG00000040229-Gpr34
Negative: ENSMUSG00000003657-Calb2, ENSMUSG00000020135-Apc2, ENSMUSG00000068740-Celsr2, ENSMUSG00000048078-Tenm4, ENSMUSG00000042436-Mfap4, ENSMUSG00000047787-Flrt1, ENSMUSG00000032181-Scg3, ENSMUSG00000027309-4930402H24Rik, ENSMUSG00000025400-Tac2, ENSMUSG00000032086-Bace1
ENSMUSG00000031517-Gpm6a, ENSMUSG00000037984-Neurod6, ENSMUSG00000023927-Satb1, ENSMUSG00000025203-Scd2, ENSMUSG00000027282-Mtch2, ENSMUSG00000022619-Mapk8ip2, ENSMUSG00000046573-Lyrm4, ENSMUSG00000027692-Tnik, ENSMUSG00000058254-Tspan7, ENSMUSG00000019124-Scrn1
ENSMUSG00000024754-Tmem2, ENSMUSG00000031617-Tmem184c, ENSMUSG00000007021-Syngr3, ENSMUSG00000040111-Gramd1b, ENSMUSG00000006333-Rps9, ENSMUSG00000031821-Gins2, ENSMUSG00000031822-Gse1, ENSMUSG00000043733-Ptpn11, ENSMUSG00000033998-Kcnk1, ENSMUSG00000000202-Btbd17
PC_ 3
Positive: ENSMUSG00000020423-Btg2, ENSMUSG00000052684-Jun, ENSMUSG00000006333-Rps9, ENSMUSG00000024190-Dusp1, ENSMUSG00000052837-Junb, ENSMUSG00000030737-Slco2b1, ENSMUSG00000031821-Gins2, ENSMUSG00000028976-Slc2a5, ENSMUSG00000036931-Nfkbid, ENSMUSG00000044786-Zfp36
ENSMUSG00000027399-Il1a, ENSMUSG00000031822-Gse1, ENSMUSG00000035356-Nfkbiz, ENSMUSG00000024073-Birc6, ENSMUSG00000026923-Notch1, ENSMUSG00000068740-Celsr2, ENSMUSG00000026628-Atf3, ENSMUSG00000000078-Klf6, ENSMUSG00000048078-Tenm4, ENSMUSG00000038393-Txnip
ENSMUSG00000019039-Dalrd3, ENSMUSG00000000028-Cdc45, ENSMUSG00000031613-Hpgd, ENSMUSG00000022817-Itgb5, ENSMUSG00000045730-Adrb2, ENSMUSG00000026965-Anapc2, ENSMUSG00000024401-Tnf, ENSMUSG00000002496-Tsc2, ENSMUSG00000027309-4930402H24Rik, ENSMUSG00000024975-Pdcd4
Negative: ENSMUSG00000092511-Gm20547, ENSMUSG00000090231-Cfb, ENSMUSG00000025498-Irf7, ENSMUSG00000074896-Ifit3, ENSMUSG00000017493-Igfbp4, ENSMUSG00000030560-Ctsc, ENSMUSG00000075014-Gm10800, ENSMUSG00000062488-Ifit3b, ENSMUSG00000035493-Tgfbi, ENSMUSG00000035385-Ccl2
ENSMUSG00000022037-Clu, ENSMUSG00000063903-Klk1, ENSMUSG00000033880-Lgals3bp, ENSMUSG00000060802-B2m, ENSMUSG00000038642-Ctss, ENSMUSG00000064373-Selenop, ENSMUSG00000033208-S100b, ENSMUSG00000002289-Angptl4, ENSMUSG00000090698-Apold1, ENSMUSG00000061232-H2-K1
ENSMUSG00000042659-Arrdc4, ENSMUSG00000036908-Unc93b1, ENSMUSG00000022257-Laptm4b, ENSMUSG00000071637-Cebpd, ENSMUSG00000036256-Igfbp7, ENSMUSG00000017390-Aldoc, ENSMUSG00000027199-Gatm, ENSMUSG00000036905-C1qb, ENSMUSG00000058715-Fcer1g, ENSMUSG00000023367-Tmem176a
PC_ 4
Positive: ENSMUSG00000035686-Thrsp, ENSMUSG00000027004-Frzb, ENSMUSG00000027878-Notch2, ENSMUSG00000038668-Lpar1, ENSMUSG00000020932-Gfap, ENSMUSG00000037206-Islr, ENSMUSG00000052957-Gas1, ENSMUSG00000021250-Fos, ENSMUSG00000029413-Naaa, ENSMUSG00000015090-Ptgds
ENSMUSG00000059325-Hopx, ENSMUSG00000028195-Cyr61, ENSMUSG00000028214-Gem, ENSMUSG00000001025-S100a6, ENSMUSG00000067818-Myl9, ENSMUSG00000052684-Jun, ENSMUSG00000021702-Thbs4, ENSMUSG00000053398-Phgdh, ENSMUSG00000028583-Pdpn, ENSMUSG00000037031-Tspan15
ENSMUSG00000032348-Gsta4, ENSMUSG00000071637-Cebpd, ENSMUSG00000024190-Dusp1, ENSMUSG00000053519-Kcnip1, ENSMUSG00000028691-Prdx1, ENSMUSG00000059970-Hspa2, ENSMUSG00000023034-Nr4a1, ENSMUSG00000021186-Fbln5, ENSMUSG00000027253-Lrp4, ENSMUSG00000020423-Btg2
Negative: ENSMUSG00000033998-Kcnk1, ENSMUSG00000060961-Slc4a4, ENSMUSG00000038094-Atp13a4, ENSMUSG00000036949-Slc39a12, ENSMUSG00000022122-Ednrb, ENSMUSG00000055003-Lrtm2, ENSMUSG00000027200-Sema6d, ENSMUSG00000003657-Calb2, ENSMUSG00000040055-Gjb6, ENSMUSG00000020135-Apc2
ENSMUSG00000037348-Paqr7, ENSMUSG00000030495-Slc7a10, ENSMUSG00000026657-Frmd4a, ENSMUSG00000053004-Hrh1, ENSMUSG00000073424-Cyp4f15, ENSMUSG00000033849-B3galt2, ENSMUSG00000001260-Gabrg1, ENSMUSG00000037697-Ddhd1, ENSMUSG00000020734-Grin2c, ENSMUSG00000035376-Hacd2
ENSMUSG00000024411-Aqp4, ENSMUSG00000096054-Syne1, ENSMUSG00000002475-Abhd3, ENSMUSG00000022537-Tmem44, ENSMUSG00000010064-Slc38a3, ENSMUSG00000047787-Flrt1, ENSMUSG00000052387-Trpm3, ENSMUSG00000028766-Alpl, ENSMUSG00000000202-Btbd17, ENSMUSG00000032036-Kirrel3
PC_ 5
Positive: ENSMUSG00000033720-Sfxn5, ENSMUSG00000058254-Tspan7, ENSMUSG00000017390-Aldoc, ENSMUSG00000005089-Slc1a2, ENSMUSG00000053398-Phgdh, ENSMUSG00000059325-Hopx, ENSMUSG00000020591-Ntsr2, ENSMUSG00000035805-Mlc1, ENSMUSG00000041556-Fbxo2, ENSMUSG00000032281-Acsbg1
ENSMUSG00000004902-Slc25a18, ENSMUSG00000028976-Slc2a5, ENSMUSG00000035686-Thrsp, ENSMUSG00000020932-Gfap, ENSMUSG00000089675-Ugt1a8, ENSMUSG00000090175-Ugt1a9, ENSMUSG00000089943-Ugt1a5, ENSMUSG00000090171-Ugt1a2, ENSMUSG00000090165-Ugt1a10, ENSMUSG00000089960-Ugt1a1
ENSMUSG00000090145-Ugt1a6b, ENSMUSG00000021792-Fam213a, ENSMUSG00000054545-Ugt1a6a, ENSMUSG00000026628-Atf3, ENSMUSG00000030737-Slco2b1, ENSMUSG00000026701-Prdx6, ENSMUSG00000027399-Il1a, ENSMUSG00000004558-Ndrg2, ENSMUSG00000045608-Dbx2, ENSMUSG00000030701-Plekhb1
Negative: ENSMUSG00000029718-Pcolce, ENSMUSG00000017344-Vtn, ENSMUSG00000025780-Itih5, ENSMUSG00000061353-Cxcl12, ENSMUSG00000042116-Vwa1, ENSMUSG00000029661-Col1a2, ENSMUSG00000036256-Igfbp7, ENSMUSG00000030717-Nupr1, ENSMUSG00000004665-Cnn2, ENSMUSG00000022548-Apod
ENSMUSG00000003617-Cp, ENSMUSG00000050212-Eva1b, ENSMUSG00000024620-Pdgfrb, ENSMUSG00000031538-Plat, ENSMUSG00000004951-Hspb1, ENSMUSG00000036814-Slc6a20a, ENSMUSG00000022090-Pdlim2, ENSMUSG00000031503-Col4a2, ENSMUSG00000027800-Tm4sf1, ENSMUSG00000063796-Slc22a8
ENSMUSG00000008999-Bmp7, ENSMUSG00000072941-Sod3, ENSMUSG00000017754-Pltp, ENSMUSG00000020467-Efemp1, ENSMUSG00000062980-Cped1, ENSMUSG00000053279-Aldh1a1, ENSMUSG00000024168-Tmem204, ENSMUSG00000001506-Col1a1, ENSMUSG00000089945-Pakap, ENSMUSG00000038729-Akap2
merged <- FindNeighbors(object = merged)
merged <- FindClusters(object = merged)
Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck
Number of nodes: 1862
Number of edges: 56945
Running Louvain algorithm...
Maximum modularity in 10 random starts: 0.8731
Number of communities: 14
Elapsed time: 0 seconds
merged <- RunTSNE(object = merged)
merged <- RunUMAP(object = merged, reduction = "pca", dims = RANGE_PCS)
p <- DimPlot(object = merged, reduction = "tsne", group.by = "plate")
print(p + theme(legend.position = "bottom"))
p2 <- DimPlot(object = merged, reduction = "umap", group.by = "plate",
label = TRUE, repel = TRUE) + NoLegend()
print(p2)
rm(merged)
Second, CCA+MNN Seurat v3 integration
for (i in 1:length(so)) {
so[[i]] <- NormalizeData(object = so[[i]], verbose = FALSE)
so[[i]] <- FindVariableFeatures(object = so[[i]], selection.method = "vst",
nfeatures = 2000, verbose = FALSE)
}
jaeger.anchors <- FindIntegrationAnchors(object.list = so, dims = RANGE_PCS)
integrated <- IntegrateData(anchorset = jaeger.anchors, dims = RANGE_PCS)
DefaultAssay(object = integrated) <- "integrated"
## Run the standard workflow for visualization and clustering
integrated <- ScaleData(object = integrated, verbose = FALSE,
display.progress = FALSE, vars.to.regress = c("nUMI", "pct_counts_mito"))
integrated <- RunPCA(object = integrated, npcs = max(RANGE_PCS),
verbose = FALSE)
integrated <- RunTSNE(object = integrated, npcs = max(RANGE_PCS),
verbose = FALSE)
ElbowPlot(integrated)
RANGE_PCS_RED <- 1:10
integrated <- RunUMAP(object = integrated, reduction = "pca",
dims = RANGE_PCS_RED)
p1 <- DimPlot(object = integrated, reduction = "tsne", group.by = "plate")
p2 <- DimPlot(object = integrated, reduction = "umap", group.by = "plate",
label = TRUE, repel = TRUE) + NoLegend()
p3 <- DimPlot(object = integrated, reduction = "pca", group.by = "plate")
plot(p3)
plot(p1)
plot(p2)
RES <- "integrated_snn_res.1"
res <- 1
integrated <- FindNeighbors(object = integrated, dims = RANGE_PCS_RED)
integrated <- FindClusters(object = integrated, resolution = res)
Modularity Optimizer version 1.3.0 by Ludo Waltman and Nees Jan van Eck
Number of nodes: 1862
Number of edges: 65447
Running Louvain algorithm...
Maximum modularity in 10 random starts: 0.7634
Number of communities: 11
Elapsed time: 0 seconds
integrated@active.assay <- "RNA"
Idents(object = integrated) <- RES
print(DimPlot(object = integrated, reduction = "umap", group.by = "plate"))
print(FeaturePlot(object = integrated, reduction = "umap", features = "pct_counts_mito"))
print(FeaturePlot(object = integrated, reduction = "umap", features = "nFeature_RNA"))
print(FeaturePlot(object = integrated, reduction = "umap", features = "nCount_RNA"))
print(DimPlot(object = integrated, reduction = "umap", group.by = RES))
markers <- FindAllMarkers(object = integrated, only.pos = TRUE,
min.pct = 0.25, thresh.use = 0.25)
Calculating cluster 0
Calculating cluster 1
Calculating cluster 2
Calculating cluster 3
Calculating cluster 4
Calculating cluster 5
Calculating cluster 6
Calculating cluster 7
Calculating cluster 8
Calculating cluster 9
Calculating cluster 10
DT::datatable(markers %>% as.data.frame() %>% dplyr::mutate_if(is.numeric,
funs(round(., 2))), extensions = c("Buttons", "FixedColumns"),
rownames = FALSE, options = list(dom = "Bfrtip", scrollX = TRUE,
fixedColumns = list(leftColumns = 1), buttons = c("csv",
"excel")))
markers_head <- as.data.frame(markers %>% group_by(cluster) %>%
top_n(n = 5, wt = avg_logFC))$gene
## graphics.off()
pheatmap(GetAssayData(object = integrated)[markers_head, ], color = viridis(100),
cluster_rows = FALSE, cluster_cols = TRUE, show_rownames = TRUE,
show_colnames = FALSE, clustering_distance_cols = "euclidean",
clustering_method = "ward.D2", clustering_distance_rows = "euclidean",
fontsize_row = 12, annotation_col = integrated@meta.data[,
c("plate", "nCount_RNA", "nFeature_RNA", "pct_counts_mito",
RES), drop = FALSE] %>% as.data.frame(), scale = "none")
save(sce, file = sprintf("%s_sce.RData", analysis_tag))
save(so, file = sprintf("%s_so.RData", analysis_tag))
date()
[1] "Thu Jul 18 14:28:30 2019"
devtools::session_info()
─ Session info ──────────────────────────────────────────────────────────
setting value
version R version 3.6.0 (2019-04-26)
os Ubuntu 16.04.5 LTS
system x86_64, linux-gnu
ui X11
language en_CA:en
collate en_CA.UTF-8
ctype en_CA.UTF-8
tz Europe/Zurich
date 2019-07-18
─ Packages ──────────────────────────────────────────────────────────────
package * version date lib source
AnnotationDbi 1.46.0 2019-05-02 [1] Bioconductor
ape 5.3 2019-03-17 [1] CRAN (R 3.6.0)
assertthat 0.2.1 2019-03-21 [1] CRAN (R 3.6.0)
backports 1.1.4 2019-04-10 [1] CRAN (R 3.6.0)
beeswarm 0.2.3 2016-04-25 [1] CRAN (R 3.6.0)
bibtex 0.4.2 2017-06-30 [1] CRAN (R 3.6.0)
Biobase * 2.44.0 2019-05-02 [1] Bioconductor
BiocGenerics * 0.30.0 2019-05-02 [1] Bioconductor
BiocNeighbors 1.2.0 2019-05-02 [1] Bioconductor
BiocParallel * 1.18.0 2019-05-03 [1] Bioconductor
BiocSingular 1.0.0 2019-05-02 [1] Bioconductor
biomaRt * 2.40.0 2019-05-02 [1] Bioconductor
bit 1.1-14 2018-05-29 [1] CRAN (R 3.6.0)
bit64 0.9-7 2017-05-08 [1] CRAN (R 3.6.0)
bitops 1.0-6 2013-08-17 [1] CRAN (R 3.6.0)
blob 1.1.1 2018-03-25 [1] CRAN (R 3.6.0)
Cairo * 1.5-10 2019-03-28 [1] CRAN (R 3.6.0)
callr 3.2.0 2019-03-15 [1] CRAN (R 3.6.0)
caTools 1.17.1.2 2019-03-06 [1] CRAN (R 3.6.0)
cli 1.1.0 2019-03-19 [1] CRAN (R 3.6.0)
cluster 2.0.9 2019-05-01 [1] CRAN (R 3.6.0)
codetools 0.2-16 2018-12-24 [3] CRAN (R 3.6.0)
colorspace 1.4-1 2019-03-18 [1] CRAN (R 3.6.0)
colourpicker 1.0 2017-09-27 [1] CRAN (R 3.6.0)
cowplot 0.9.4 2019-01-08 [1] CRAN (R 3.6.0)
crayon 1.3.4 2017-09-16 [1] CRAN (R 3.6.0)
crosstalk 1.0.0 2016-12-21 [1] CRAN (R 3.6.0)
data.table 1.12.2 2019-04-07 [1] CRAN (R 3.6.0)
DBI 1.0.0 2018-05-02 [1] CRAN (R 3.6.0)
DelayedArray * 0.10.0 2019-05-02 [1] Bioconductor
DelayedMatrixStats 1.6.0 2019-05-02 [1] Bioconductor
desc 1.2.0 2018-05-01 [1] CRAN (R 3.6.0)
devtools 2.0.2 2019-04-08 [1] CRAN (R 3.6.0)
digest 0.6.19 2019-05-20 [1] CRAN (R 3.6.0)
dplyr * 0.8.1 2019-05-14 [1] CRAN (R 3.6.0)
dqrng 0.2.1 2019-05-17 [1] CRAN (R 3.6.0)
DT * 0.5 2018-11-05 [1] CRAN (R 3.6.0)
dynamicTreeCut 1.63-1 2016-03-11 [1] CRAN (R 3.6.0)
edgeR 3.26.5 2019-06-21 [1] Bioconductor
evaluate 0.14 2019-05-28 [1] CRAN (R 3.6.0)
fitdistrplus 1.0-14 2019-01-23 [1] CRAN (R 3.6.0)
formatR 1.7 2019-06-11 [1] CRAN (R 3.6.0)
fs 1.3.0 2019-05-02 [1] CRAN (R 3.6.0)
future 1.12.0 2019-03-08 [1] CRAN (R 3.6.0)
future.apply 1.2.0 2019-03-07 [1] CRAN (R 3.6.0)
gbRd 0.4-11 2012-10-01 [1] CRAN (R 3.6.0)
gdata 2.18.0 2017-06-06 [1] CRAN (R 3.6.0)
GenomeInfoDb * 1.20.0 2019-05-02 [1] Bioconductor
GenomeInfoDbData 1.2.1 2019-05-04 [1] Bioconductor
GenomicRanges * 1.36.0 2019-05-02 [1] Bioconductor
ggbeeswarm 0.6.0 2017-08-07 [1] CRAN (R 3.6.0)
ggplot2 * 3.2.0 2019-06-16 [1] CRAN (R 3.6.0)
ggrepel 0.8.0 2018-05-09 [1] CRAN (R 3.6.0)
ggridges 0.5.1 2018-09-27 [1] CRAN (R 3.6.0)
globals 0.12.4 2018-10-11 [1] CRAN (R 3.6.0)
glue 1.3.1 2019-03-12 [1] CRAN (R 3.6.0)
gplots 3.0.1.1 2019-01-27 [1] CRAN (R 3.6.0)
gridExtra 2.3 2017-09-09 [1] CRAN (R 3.6.0)
gtable 0.3.0 2019-03-25 [1] CRAN (R 3.6.0)
gtools 3.8.1 2018-06-26 [1] CRAN (R 3.6.0)
hms 0.4.2 2018-03-10 [1] CRAN (R 3.6.0)
htmltools 0.3.6 2017-04-28 [1] CRAN (R 3.6.0)
htmlwidgets 1.3 2018-09-30 [1] CRAN (R 3.6.0)
httpuv 1.5.1 2019-04-05 [1] CRAN (R 3.6.0)
httr 1.4.0 2018-12-11 [1] CRAN (R 3.6.0)
ica 1.0-2 2018-05-24 [1] CRAN (R 3.6.0)
igraph 1.2.4.1 2019-04-22 [1] CRAN (R 3.6.0)
IRanges * 2.18.1 2019-05-31 [1] Bioconductor
irlba 2.3.3 2019-02-05 [1] CRAN (R 3.6.0)
iSEE * 1.4.0 2019-05-02 [1] Bioconductor
jsonlite 1.6 2018-12-07 [1] CRAN (R 3.6.0)
KernSmooth 2.23-15 2015-06-29 [3] CRAN (R 3.6.0)
knitr * 1.23 2019-05-18 [1] CRAN (R 3.6.0)
labeling 0.3 2014-08-23 [1] CRAN (R 3.6.0)
later 0.8.0 2019-02-11 [1] CRAN (R 3.6.0)
lattice 0.20-38 2018-11-04 [3] CRAN (R 3.6.0)
lazyeval 0.2.2 2019-03-15 [1] CRAN (R 3.6.0)
limma 3.40.2 2019-05-17 [1] Bioconductor
listenv 0.7.0 2018-01-21 [1] CRAN (R 3.6.0)
lmtest 0.9-37 2019-04-30 [1] CRAN (R 3.6.0)
locfit 1.5-9.1 2013-04-20 [1] CRAN (R 3.6.0)
lsei 1.2-0 2017-10-23 [1] CRAN (R 3.6.0)
magrittr 1.5 2014-11-22 [1] CRAN (R 3.6.0)
MASS 7.3-51.4 2019-03-31 [3] CRAN (R 3.6.0)
Matrix * 1.2-17 2019-03-22 [1] CRAN (R 3.6.0)
matrixStats * 0.54.0 2018-07-23 [1] CRAN (R 3.6.0)
memoise 1.1.0 2017-04-21 [1] CRAN (R 3.6.0)
metap 1.1 2019-02-06 [1] CRAN (R 3.6.0)
mgcv 1.8-28 2019-03-21 [3] CRAN (R 3.6.0)
mime 0.7 2019-06-11 [1] CRAN (R 3.6.0)
miniUI 0.1.1.1 2018-05-18 [1] CRAN (R 3.6.0)
munsell 0.5.0 2018-06-12 [1] CRAN (R 3.6.0)
nlme 3.1-139 2019-04-09 [3] CRAN (R 3.6.0)
npsurv 0.4-0 2017-10-14 [1] CRAN (R 3.6.0)
pbapply 1.4-0 2019-02-05 [1] CRAN (R 3.6.0)
pheatmap * 1.0.12 2019-01-04 [1] CRAN (R 3.6.0)
pillar 1.4.1 2019-05-28 [1] CRAN (R 3.6.0)
pkgbuild 1.0.3 2019-03-20 [1] CRAN (R 3.6.0)
pkgconfig 2.0.2 2018-08-16 [1] CRAN (R 3.6.0)
pkgload 1.0.2 2018-10-29 [1] CRAN (R 3.6.0)
plotly 4.9.0 2019-04-10 [1] CRAN (R 3.6.0)
plyr 1.8.4 2016-06-08 [1] CRAN (R 3.6.0)
png 0.1-7 2013-12-03 [1] CRAN (R 3.6.0)
prettyunits 1.0.2 2015-07-13 [1] CRAN (R 3.6.0)
processx 3.3.0 2019-03-10 [1] CRAN (R 3.6.0)
progress 1.2.2 2019-05-16 [1] CRAN (R 3.6.0)
promises 1.0.1 2018-04-13 [1] CRAN (R 3.6.0)
ps 1.3.0 2018-12-21 [1] CRAN (R 3.6.0)
purrr 0.3.2 2019-03-15 [1] CRAN (R 3.6.0)
R.methodsS3 1.7.1 2016-02-16 [1] CRAN (R 3.6.0)
R.oo 1.22.0 2018-04-22 [1] CRAN (R 3.6.0)
R.utils 2.8.0 2019-02-14 [1] CRAN (R 3.6.0)
R6 2.4.0 2019-02-14 [1] CRAN (R 3.6.0)
RANN 2.6.1 2019-01-08 [1] CRAN (R 3.6.0)
RColorBrewer 1.1-2 2014-12-07 [1] CRAN (R 3.6.0)
Rcpp 1.0.1 2019-03-17 [1] CRAN (R 3.6.0)
RCurl 1.95-4.12 2019-03-04 [1] CRAN (R 3.6.0)
Rdpack 0.11-0 2019-04-14 [1] CRAN (R 3.6.0)
remotes 2.0.4 2019-04-10 [1] CRAN (R 3.6.0)
rentrez 1.2.2 2019-05-02 [1] CRAN (R 3.6.0)
reshape2 1.4.3 2017-12-11 [1] CRAN (R 3.6.0)
reticulate 1.12 2019-04-12 [1] CRAN (R 3.6.0)
rintrojs 0.2.0 2017-07-04 [1] CRAN (R 3.6.0)
rlang 0.4.0 2019-06-25 [1] CRAN (R 3.6.0)
rmarkdown 1.12 2019-03-14 [1] CRAN (R 3.6.0)
ROCR 1.0-7 2015-03-26 [1] CRAN (R 3.6.0)
rprojroot 1.3-2 2018-01-03 [1] CRAN (R 3.6.0)
RSQLite 2.1.1 2018-05-06 [1] CRAN (R 3.6.0)
rsvd 1.0.1 2019-06-02 [1] CRAN (R 3.6.0)
Rtsne 0.15 2018-11-10 [1] CRAN (R 3.6.0)
S4Vectors * 0.22.0 2019-05-02 [1] Bioconductor
scales 1.0.0 2018-08-09 [1] CRAN (R 3.6.0)
scater * 1.12.2 2019-05-24 [1] Bioconductor
scran * 1.12.1 2019-05-27 [1] Bioconductor
sctransform 0.2.0 2019-04-12 [1] CRAN (R 3.6.0)
SDMTools 1.1-221.1 2019-04-18 [1] CRAN (R 3.6.0)
sessioninfo 1.1.1 2018-11-05 [1] CRAN (R 3.6.0)
Seurat * 3.0.0 2019-04-15 [1] CRAN (R 3.6.0)
shiny 1.3.2 2019-04-22 [1] CRAN (R 3.6.0)
shinyAce 0.3.3 2019-01-03 [1] CRAN (R 3.6.0)
shinydashboard 0.7.1 2018-10-17 [1] CRAN (R 3.6.0)
shinyjs 1.0 2018-01-08 [1] CRAN (R 3.6.0)
SingleCellExperiment * 1.6.0 2019-05-02 [1] Bioconductor
statmod 1.4.32 2019-05-29 [1] CRAN (R 3.6.0)
stringi 1.4.3 2019-03-12 [1] CRAN (R 3.6.0)
stringr 1.4.0 2019-02-10 [1] CRAN (R 3.6.0)
SummarizedExperiment * 1.14.0 2019-05-02 [1] Bioconductor
survival 2.44-1.1 2019-04-01 [3] CRAN (R 3.6.0)
testthat 2.1.1 2019-04-23 [1] CRAN (R 3.6.0)
tibble 2.1.3 2019-06-06 [1] CRAN (R 3.6.0)
tidyr 0.8.3 2019-03-01 [1] CRAN (R 3.6.0)
tidyselect 0.2.5 2018-10-11 [1] CRAN (R 3.6.0)
tsne 0.1-3 2016-07-15 [1] CRAN (R 3.6.0)
usethis 1.5.0 2019-04-07 [1] CRAN (R 3.6.0)
vipor 0.4.5 2017-03-22 [1] CRAN (R 3.6.0)
viridis * 0.5.1 2018-03-29 [1] CRAN (R 3.6.0)
viridisLite * 0.3.0 2018-02-01 [1] CRAN (R 3.6.0)
withr 2.1.2 2018-03-15 [1] CRAN (R 3.6.0)
xfun 0.8 2019-06-25 [1] CRAN (R 3.6.0)
XML 3.98-1.20 2019-06-06 [1] CRAN (R 3.6.0)
xtable 1.8-4 2019-04-21 [1] CRAN (R 3.6.0)
XVector 0.24.0 2019-05-02 [1] Bioconductor
yaml 2.2.0 2018-07-25 [1] CRAN (R 3.6.0)
zlibbioc 1.30.0 2019-05-02 [1] Bioconductor
zoo 1.8-5 2019-03-21 [1] CRAN (R 3.6.0)
[1] /home/Shared/Rlib/release-3.9-lib
[2] /home/imallona/R/x86_64-pc-linux-gnu-library/3.6
[3] /usr/local/R/R-3.6.0/library